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Abstract 

We generalize lattice models of brittle fracture to arbitrary nonlinear force 
laws and study the existence of arrested semi-infinite cracks. Unlike what is 
seen in the discontinuous case studied to date, the range in driving displace- 
ment for which these arrested cracks exist is either very small or precisely 
zero. Also, our results indicate that small changes in the vicinity of the crack 
tip can have an extremely large effect on arrested cracks. Finally, we briefly 
discuss the possible relevance of our findings to recent experiments. 
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Recent years have seen a rebirth of interest by the physics community in the issue of 
dynamic fracture. This is due to a variety of new experimental results which are not explain- 
able within the confines of the traditional engineering approach to fracture These results 
include a dynamical instability to micro-branching the formation of non-smooth frac- 
ture surfaces ||] and the rapid variation of the fracture energy (including dissipative losses 
incurred during cleavage) with crack velocity [gj. These issues are reviewed in a recent paper 
by Fineberg and Marder ||. 

One approach for dealing with dynamic fracture involves restricting the atomic interac- 
tions to those occurring between neighboring sites of an originally unstrained lattice. These 
lattice models can never be as realistic as full molecular dynamics simulations, but compen- 
sate for this shortcoming by being much more amenable to analysis, both numerical and 
(via the Wiener-Hopf technique) otherwise. This approach was pioneered by Slepyan and 
co-workers and further developed by Marder and Gross || and most recently by our- 
selves ||. Most of the results to date have been obtained using a simplified force law which 
is linear until some threshold displacement at which point it drops abruptly to zero. Below, 
we will study a generalization for which the force is a smooth function of the lattice strain. 
One of our goals is to learn which aspects of fracture are sensitive to microscopic details and 
which are universal. 

One interesting aspect of these lattice models concern the existence of a range of driving 
displacements A for which non-moving semi-infinite crack solutions can be found. For 
the aforementioned discontinuous force model, there exists a wide range of these arrested 
cracks. For example, ref. |9j found that A could range from 40% below to 40% above 
the Griffith displacement Ac, the driving at which it first becomes energetically favorable 
for the system to crack. This phenomena is connected to the existence of a velocity gap, 
i.e. a minimal velocity for stable crack propagation. Experimentally, no such gap has 
been reported, even for materials such as single- crystal silicon [T(| which should be at least 
approximately describable by lattice models. It is therefore of some interest to study how 
the arrested crack range depends on the microscopic details of the assumed atomic force 
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law. Here we present the results of such as a study, including the finding that this range 
drops rapidly towards zero as the force law is made smoother and hence more realistic. 

As in ref. ||, we work with a square lattice and with scalar displacements (mode III). 
We focus on arrested cracks and write the static equation as 

= -/ (U i+ ij — Uij) + f (Uij —Ui-tj) - f (Uij + i-Uij) + f (mj-Uij-i) (1) 

Here the indices {i,j} label the lattice site and u is the displacement. Sites on the last row 
of the the lattice, j = N y , are coupled to a row with fixed displacement A. The first row, 
j = 1, is coupled to a j = displacement field u ij0 which via symmetry equals — u^i. Finally, 
/ is a nonlinear function of its argument, the lattice strain. We investigate two forms [|TT|1 : 

1 + tanh(q(l - u)) 

JeW = ~U 2 

1 + tanh a 

ua a+1 

fju) = — - r— —r (3) 

For both of these forms, increasing a reduces the length scale over which / falls to zero once 
outside the Hooke's law regime (u < 1). The exponential force f e reduces to the familiar 
discontinuous force (linear until complete failure) as a — > oo. 

Our procedure for finding solutions is in principle straightforward. At large positive i 
in the uncracked material, we know that the system will adopt a uniformly strained state. 
Conversely, at large negative i the cracked state will have a large displacement it^i and 
(almost) zero strains for j > 1. Fixing the boundary condition A allows us to easily find 
these asymptotic states. Once found, these solutions are used as fixed displacements for the 
columns i = N x + 1 and i = —N x — 1 respectively. The arrested crack then requires us to 
solve for (2N X + l)N y variables. We impose the equation at motion at all sites except for the 
crack "tip" , (i — 0, j — 1) where instead we specify the displacement; this approach perserves 
the banded structure of the system. Newton's algorithm then allows us to converge to a 
solution. Afterwards, the residual equation of motion becomes a solvability condition with 
which A can be determined. The range of allowed values of A for arrested cracks is found 
as one systematically sweeps through the value of the aforementioned fixed displacement. 



In Fig. 1 we present our results for the exponential model. For illustration, we have 
chosen to show data for N y = 10 as a function of a. For large a, the range of A is large and 
there is a marked asymmetry between the rising segment of A versus imposed displacement 
and the (much steeper) falling segment. As a — > oo, the falling portion becomes vertical. 
These segments represent different crack solutions at fixed A; as A reaches the end of 
its allowed range, these solution branches collide and disappear in a standard saddle-node 
bifurcation point. To verify this, we have performed |12[] a linear stability calculation of 
these solutions, assuming purely inertial dynamics (i.e. setting the left hand side of Eq. 1 to 
Uij). As expected, there is a single mode of the spectrum for the growth rate to for which uj 2 
goes from negative to positive as we go up the rising segment, reach the maximal driving, 
and then go back down. 
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FIG. 1. A versus imposed Uq i displacement for different values of a in f e . All data has N y = 10 
and N x = 100. 

Fig. 1 demonstrates that as the potential is made smoother, the range of arrested cracks 
shrinks dramatically. In Fig. 2, we show this range as a percentage of Ac- The best fit to 
our data suggests that the range vanishes as an essentially singular function of a, 
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where for N y = 10, a ~ 6.6 and otherwise is a slowly varying function of N y as long as the 
system is sufficiently large compared to the potential fall-off. 
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FIG. 2. Arrested crack range normalized by the Griffith displacement Aq versus a in f e ; again 
all data is for N y = 10, N x = 100. 

Let us now turn to the power-law form. Based on our findings above, we would expect 
that this rather smooth force law would give rise to a range which is practically zero. We 
have verified this prediction in two ways. First, for the case a = 3 we performed our usual 
scan over imposed w ,i displacement and noted that the selected A varies by less than 10~ 6 . 
Second, we computed the stability spectrum and found a mode at u 2 < 10~ 6 ; this value is 
indicative of how close we are at a randomly chosen displacement to the extremal value of A 
at the saddle-node bifurcation . These numbers are consistent with our numerical accuracy 
and hence the true range is probably even smaller. Needless to say, ranges of this size would 
be unmeasurable. It is interesting to point out that the almost-zero mode is nothing other 
than a spatial translation of the crack. That is, translating the crack with respect to the 
underlying fixed lattice is almost a symmetry of the solution. 



So, by making the potential smoother one tends to eliminate arrested crack solutions. 
How does this change come about? To try to address this question, we plot in Fig. 3 
the lattice strain field iti+i,i — u^x for — N x < % > N x for the three potentials, exponential 
with a = 5 or 2 and power-law with a = 3. For this comparison, we have found (stable) 
solutions with w ,i = .75 for all three potentials, and then normalized the strains by dividing 
with the respective values of A. First, we note that beyond x ~ 5, the different cases 
are virtually indistinguishable and all lie on the expected x _1//2 universal curve [|]J§. The 
interior "process-zone" region is affected by changing the potential, but rather minimally. 
For example, the two exponential cases differ in only one or two points, yet this is sufficient 
to shrink the arrested crack range by almost an order of magnitude. The power-law choice 
has a process-zone which is a bit wider and there is less maximal strain, but that is all. We 
thus conclude that the existence and size of the arrested crack range are extremely sensitive 
to microscopic details! We note in passing that the process-zone for any specific potential 
quickly reaches an asymptotic size once N y is sufficiently large and in particular does not 
increase indefinitely in the macroscopic limit. Treatments [0,0 which include a mesoscopic 
size "cohesive-zone" are therefore not accurate representations of this class of lattice models. 



In a recent experiment 1X0 1 on fracture in silicon, no arrested cracks were observed. A 
molecular dynamics simulation using a modified Stillinger- Weber potential also exhibited 
no arrested cracks when studied at high enough temperature. However, the potentials used 
here were rather short-ranged, as compared with some estimates that arise from density- 
functional theory [|TjJ. Our results indicate that increasing the range and thereby using 
smoother potentials will eliminate (at least as far as experimentally attainable precision is 
occurred) arrested cracks and may offer a simpler explanation of the experimental finding 
than one which requires thermal creep. This could of course be tested in principle by re-doing 
the experiments at a reduced temperature. 
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FIG. 3. Strain Uj+i ; i — Ui t i for three different potentials. Data is for N y = 40, N x = 200 and is 
normalized to the respective A values. The power-law curve has been shifted two sites to the left 
so as to better match the field at large positive i. 
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